Academic Citation

If you use this code in your work or research, we kindly request that you cite our publication:

Xiaofan Lu, et al. (2025). FigureYa: A Standardized Visualization Framework for Enhancing Biomedical Data Interpretation and Research Efficiency. iMetaMed. https://doi.org/10.1002/imm3.70005

需求描述

实现利用USCS上的pan-cancer数据分别通过EPIC和Cibersort算法，计算各个肿瘤的免疫细胞浸润信息，并输出结果，同时计算某个指定基因在各个肿瘤中分别与免疫细胞的相关性，并分别输出相关性热图和相关性信息的excel表

例文只使用了一种方法，我想用多种方法来证明相关性

Requirement Description

Use the pan-cancer data on USCS to calculate the immune cell infiltration information of each tumor through the EPIC and Cibersort algorithms, and output the results, and at the same time calculate the correlation between a specified gene and immune cells in each tumor, and output the correlation heat map and correlation information in the excel sheet.

The example uses only one method, and I want to use multiple methods to prove correlation

出自https://www.frontiersin.org/articles/10.3389/fonc.2021.634617/full fromhttps://www.frontiersin.org/articles/10.3389/fonc.2021.634617/full

FIGURE 7 | Associations of CD96 expression to tumor purity and immune infiltration. (B) The correlations of CD96 expression and immune infiltration in cancers.

Besides, we employed TIMER2.0 to exhibit the landscape of CD96 correlating with various immune infiltrates in human cancers (Figure 7B). # 应用场景

分别用TIMER、EPIC、Cibersort计算各个肿瘤的免疫细胞浸润信息。

热图可用于展示目的基因在泛癌中与免疫浸润之间的关系。

下面分别展示TIMER、EPIC、Cibersort这三种计算方法并画图（例文用的是TIMER2.0）。三种方法的细胞不同，可相互参照，结合研究兴趣，讨论自己的结果。

Application Scenarios

TIMER, EPIC, and Cibersort were used to calculate the immune cell infiltration information of each tumor.

Heatmaps can be used to demonstrate the relationship between genes of interest and immune invasion in pan-cancer.

The following three calculation methods are shown and plotted respectively (TIMER 2.0 is used in the example). The cells of the three methods are different, and they can be cross-referenced, combined with research interests, and discuss their own results.

环境设置

Environment settings

source("install_dependencies.R")

## 开始安装R包...
## ===========================================
## 
## 安装CRAN包...
## 包已安装: caret 
## 包已安装: devtools 
## 包已安装: e1071 
## 包已安装: estimate 
## 包已安装: ggplot2 
## 包已安装: ggpubr 
## 包已安装: glmnet 
## 包已安装: limSolve 
## 包已安装: patchwork 
## 包已安装: ppcor 
## 包已安装: pracma 
## 包已安装: preprocessCore 
## 包已安装: showtext 
## 包已安装: survival 
## 包已安装: survminer 
## 包已安装: sva 
## 包已安装: tibble 
## 包已安装: tidyHeatmap 
## 包已安装: timeROC 
## 
## 安装Bioconductor包...
## 包已安装: DESeq2 
## 包已安装: biomaRt 
## 包已安装: limma 
## 包已安装: EPIC 
## 
## ===========================================
## 包安装完成！
## 现在可以运行此目录中的R脚本了。

library(ggplot2)

## Warning: package 'ggplot2' was built under R version 4.4.3

library(ggpubr)
library(patchwork)
library(showtext)

## Loading required package: sysfonts

## Loading required package: showtextdb

library(EPIC)
library(IOBR)

## Loading required package: tibble

## Loading required package: dplyr

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

## Loading required package: survival

## Loading required package: ComplexHeatmap

## Loading required package: grid

## ========================================
## ComplexHeatmap version 2.20.0
## Bioconductor page: http://bioconductor.org/packages/ComplexHeatmap/
## Github page: https://github.com/jokergoo/ComplexHeatmap
## Documentation: http://jokergoo.github.io/ComplexHeatmap-reference
## 
## If you use it in published research, please cite either one:
## - Gu, Z. Complex Heatmap Visualization. iMeta 2022.
## - Gu, Z. Complex heatmaps reveal patterns and correlations in multidimensional 
##     genomic data. Bioinformatics 2016.
## 
## 
## The new InteractiveComplexHeatmap package can directly export static 
## complex heatmaps into an interactive Shiny app with zero effort. Have a try!
## 
## This message can be suppressed by:
##   suppressPackageStartupMessages(library(ComplexHeatmap))
## ========================================

## Loading required package: tidyHeatmap

## Warning: package 'tidyHeatmap' was built under R version 4.4.2

## ========================================
## tidyHeatmap version 1.8.1
## If you use tidyHeatmap in published research, please cite:
## 1) Mangiola et al. tidyHeatmap: an R package for modular heatmap production 
##   based on tidy principles. JOSS 2020.
## 2) Gu, Z. Complex heatmaps reveal patterns and correlations in multidimensional 
##   genomic data. Bioinformatics 2016.
## This message can be suppressed by:
##   suppressPackageStartupMessages(library(tidyHeatmap))
## ========================================

## 
## Attaching package: 'tidyHeatmap'

## The following object is masked from 'package:stats':
## 
##     heatmap

## Loading required package: clusterProfiler

##

## clusterProfiler v4.12.6 Learn more at https://yulab-smu.top/contribution-knowledge-mining/
## 
## Please cite:
## 
## S Xu, E Hu, Y Cai, Z Xie, X Luo, L Zhan, W Tang, Q Wang, B Liu, R Wang,
## W Xie, T Wu, L Xie, G Yu. Using clusterProfiler to characterize
## multiomics data. Nature Protocols. 2024, doi:10.1038/s41596-024-01020-z

## 
## Attaching package: 'clusterProfiler'

## The following object is masked from 'package:stats':
## 
##     filter

## Loading required package: survminer

## Warning: package 'survminer' was built under R version 4.4.2

## 
## Attaching package: 'survminer'

## The following object is masked from 'package:survival':
## 
##     myeloma

## Warning: package 'limSolve' was built under R version 4.4.2

## ==========================================================================
##   IOBR v0.99.0  Immuno-Oncology Biological Research 
##   For Tutorial: https://iobr.github.io/book/
##   For Help: https://github.com/IOBR/IOBR/issues
## 
##  If you use IOBR in published research, please cite:
##  DQ Zeng, YR Fang, WJ Qiu, ..., GC Yu*, WJ Liao*, (2024) 
##  IOBR2: Multidimensional Decoding Tumor Microenvironment for Immuno-Oncology Research. 
##  bioRxiv, 2024.01.13.575484; 
##  https://www.biorxiv.org/content/10.1101/2024.01.13.575484v2.full.pdf 
##  Higly Cited Paper and Hot Paper of WOS
## ==========================================================================

## 
## Attaching package: 'IOBR'

## The following object is masked from 'package:EPIC':
## 
##     EPIC

showtext.auto(enable = TRUE)

## 'showtext.auto()' is now renamed to 'showtext_auto()'
## The old version still works, but consider using the new function in future code

font.add('arial', 'arial.ttf') #设置字体，windows用户  # Set fonts for Windows users

## 'font.add()' is now renamed to 'font_add()'
## The old version still works, but consider using the new function in future code

Sys.setenv(LANGUAGE = "en") #显示英文报错信息 #error messages are displayed in English
options(stringsAsFactors = FALSE) #禁止chr转成factor # chr is not allowed to be converted to factor

TIMER

从TIMER下载数据

Methods: Tumor IMmune Estimation Resource 2.0 (TIMER2.0; http://timer.cistrome.org/) web server is a comprehensive resource forsystematical analysis of immune infiltrates across diverse cancertypes. At first, we used it to study the differential expressionof CD96 between tumor and adjacent normal tissues across all TCGA cohorts.

We then explored the association between CD96 expression and immune infiltration based on several immune deconvolution algorithms

进入 TIMER2.0 或者 TIMER2.0 (comp-genomics.org)，按一下步骤获得分析结果文件：

点击Immune
在Gene Expression里选择CD96
在Immune Infiltrates里每次选一种免疫细胞
点击Submit
点击Table下载分析结果
重复以上操作，直到获得所有免疫细胞的分析结果，放到data文件夹下。（虽然一次可以选择多个细胞，但网站对数量有限制，所以我们还是一次只选一种细胞。 !!!Querying 115 interested features exceeds the limited amount (n=20). Please reduce the number. ）

TIMER

Download data from TIMER

We then explored the association between CD96 expression and immune infiltration based on several immune deconvolution algorithms

Go to TIMER2.0 or TIMER2.0 (comp-genomics.org) and follow the steps below to obtain the analysis result file:

Click Immune
Select ‘CD96’ in ‘Gene Expression’
Select one type of immune cell at a time in ‘Immune Infiltrates’
Click ‘Submit’
Click ‘Table’ to download the analysis results
Repeat the above until you have the results of the analysis of all immune cells and put them in the data folder. (Although you can select more than one cell at a time, the website has a limit on the number, so we still only select one cell at a time.) !! Querying 115 interested features exceeds the limited amount (n=20). Please reduce the number. ）

批量读入TIMER结果

Batch reading TIMER results

data文件夹里除了从TIMER下载的Table表以外，不要放其他任何文件。 Do not put any files in the data folder except for the table table downloaded from the TIMER.

setwd("./data")
file <- dir()
data <- list()
for (i in (1:length(file))) {
  data[[i]] <- read.csv(file[i],header = T)
}

画图前的数据预处理

TIMER的免疫细胞分类细致，有19类，例文给每一类画成一个小热图，排成了两行，因此画图时麻烦些。

Pre-processing of data before drawing

There are 19 types of immune cells in TIMER, and each type is drawn into a small heat map arranged in two rows, so it is more troublesome to draw the picture.

# y轴排列顺序
# Y-axis sort order
name <- c("UVM (n=80)", "UCS (n=57)", "UCEC (n=545)", "THYM (n=120)", 
        "THCA (n=509)", "TGCT (n=150)", "STAD (n=415)", "SKCM (n=471)", 
        "SKCM-Primary (n=103)", "SKCM-Metastasis (n=368)", "SARC (n=260)", 
        "READ (n=166)", "PRAD (n=498)", "PCPG (n=181)", "PAAD (n=179)", 
        "OV (n=303)", "MESO (n=87)", "LUSC (n=501)", "LUAD (n=515)", 
        "LIHC (n=371)", "LGG (n=516)", "KIRP (n=290)", "KIRC (n=533)", 
        "KICH (n=66)", "HNSC (n=522)", "HNSC-HPV +  (n=98)", "HNSC-HPV- (n=422)", 
        "GBM (n=153)", "ESCA (n=185)", "DLBC (n=48)", "COAD (n=458)", 
        "CHOL (n=36)", "CESC (n=306)", "BRCA (n=1100)", "BRCA-LumB (n=219)", 
        "BRCA-LumA (n=568)", "BRCA-Her2 (n=82)", "BRCA-Basal (n=191)", 
        "BLCA (n=408)", "ACC (n=79)")

# 定义函数画heatmap
# Define the function to draw a heatmap
plot_tme <- function(x){
  x$pvalue = ifelse(x$adj.p >= 0.05, "p≥0.05", "p<0.05")
  x$cancer = factor(x$cancer, levels = name) 
  
  ggplot(x, aes(infiltrates, cancer, 
                shape = pvalue, 
                color = rho)) +
    geom_point(size = 3) + 
    scale_shape_manual(values = c(15, 7)) + 
    scale_color_gradient2(low = "#2b8cbe",mid = "white",high = "#e41a1c") + 
    #scale_x_discrete(position = "top") + 
    theme_bw() + 
    theme(axis.title.y = element_blank(),
          axis.title.x = element_blank(),
          axis.text.x = element_text(angle = 90,hjust = 0,vjust = 0),
          axis.text.x.bottom = element_text(family = "arial"))
}

# 画出19个热图
# Draw 19 heat maps
dd <- lapply(data, plot_tme)
#dd[[1]] # 可以这样查看每一个heatmap # You can view each heatmap like this

画图

Drawing

# 设置和下半部分分割的位置
# Set the position of the division with the lower half
breakpoint <- 8 # 上面画8个热图，下面画19-8=11个 # Draw 8 heat maps on the top and 19-8=11 on the bottom

## 先画上半部分
# 画第1个热图（第1个热图左侧带label，其余不带）
## Draw the top half first
# Draw the first heatmap (the left side of the first heatmap has a label, and the rest is not)
up <- dd[[1]] + 
  scale_x_discrete(position = "top") + 
  guides(color=FALSE) + guides(shape=FALSE) # 不画图例 # No legend

## Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
## of ggplot2 3.3.4.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

# 加上第2-8个热图
# Plus the 2nd-8th heatmap
for (i in 2:(breakpoint)) {
  up <- up + dd[[i]] + 
    scale_x_discrete(position = "top") + 
    theme(axis.text.y = element_blank()) + # 不画y轴label # Don't draw the y-axis label
    guides(color = FALSE) + guides(shape = FALSE) # 不画图例 # No legend
}

# 按照x轴来调整图片组合宽度
# Adjust the width of the image combination according to the x-axis
width_up <- NULL
for (i in 1:breakpoint) {
  width_up <- c(width_up,length(unique(data[[i]]$infiltrates)))
}

# 画完上半部分
# Finish the first half
up + plot_layout(guides = 'collect',widths = width_up)

## Warning: Removed 10 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 8 rows containing missing values or values outside the scale range
## (`geom_point()`).

## 画下半部分
# 画下面第1个
## Draw the bottom half
# Draw the first one below
down <- dd[[breakpoint + 1]] + scale_x_discrete(position = "bottom") + 
  guides(color = FALSE) + guides(shape = FALSE) + # 不画图例
  theme(axis.text.x = element_text(angle = 90,hjust = 1,vjust = 0))

# 画第10-18个热图
# Draw the 10th-18th heatmap
for (i in (breakpoint + 2):(length(file)-1)) {
  down <- down + dd[[i]] + scale_x_discrete(position = "bottom") + 
    theme(axis.text.y = element_blank(),
          axis.text.x = element_text(angle = 90,hjust = 1,vjust = 0)) + 
    guides(color=FALSE) + guides(shape=FALSE) # 不画图例
}

# 画第19个，画图例
# Draw the 19th one, draw a legend
down <- down + dd[[length(file)]] + scale_x_discrete(position = "bottom") + 
  theme(axis.text.y = element_blank(),
        axis.text.x = element_text(angle = 90,hjust = 1,vjust = 0))

down <- down & scale_fill_continuous(limits = c(-1, 1))

# 按照x轴来调整图片组合宽度
# Adjust the width of the image combination according to the x-axis
width_down <- NULL
for (i in (breakpoint + 1):length(file)) {
  width_down <- c(width_down,length(unique(data[[i]]$infiltrates)))
}

# 画完下半部分
# Finish the bottom half
down + plot_layout(guides = 'collect',widths = width_down)

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_point()`).

## 上下组合
## Combination of top and bottom
(up + plot_layout(guides = 'collect',widths = width_up))-(down + plot_layout(guides = 'collect',widths = width_down)) + plot_layout(nrow = 2)

## Warning: Removed 10 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 8 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_point()`).

# 保存到文件
# Save to file
ggsave("TIMER.pdf", 
       height=15, width=10)

## Warning: Removed 10 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 8 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).

## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_point()`).

# 输出的PDF文件是矢量图，可以用illustrator等矢量图编辑器打开，进一步调整图形、文字
# The output PDF file is a vector image, which can be opened with a vector editor such as illustrator to further adjust the graphics and text

EPIC

IOBR is highly recommended for immunoinfiltration algorithms （GitHub - IOBR/IOBR: IOBR is an R package to perform comprehensive analysis of tumor microenvironment and signatures for immuno-oncology.）

输入数据

pancancer_mRNA_exprSet.Rdata，表达矩阵和癌症名type

expr_pancancer.Rda，表达矩阵，取自pancancer_mRNA_exprSet.Rdata。

数据来源于UCSC xena Pan-cancer，从数据下载到获得以上两个文件的过程见TCGA_pan.R（位于当前文件夹）

Enter data

pancancer_mRNA_exprSet.Rdata, expression matrix and cancer name type

expr_pancancer. Rda, Expression Matrix, taken from pancancer_mRNA_exprSet.Rdata.

The data is from UCSC xena Pan-cancer, and the process from data download to obtaining the above two files is described in ‘TCGA_pan. R’ (located in the current folder).

# 载入数据
# Load data
(load("expr_pancancer.Rda"))

## [1] "expr"

expr[1:5,1:5]

##            TCGA-OR-A5J1 TCGA-OR-A5J2 TCGA-OR-A5J3 TCGA-OR-A5J5 TCGA-OR-A5J6
## OR4F5           -9.9658      -9.9658      -9.9658      -9.9658      -9.9658
## FO538757.3      -9.9658      -9.9658      -3.4580      -9.9658      -3.6259
## FO538757.2       4.6502       3.9783       3.5803       3.6916       4.8724
## OR4F29          -9.9658      -6.5064      -9.9658      -3.8160      -9.9658
## OR4F16          -9.9658      -6.5064      -9.9658      -3.8160      -9.9658

dim(expr)

## [1] 17146  9712

(load("pancancer_mRNA_exprSet.Rdata"))

## [1] "mRNA_exprSet"

mRNA_exprSet[1:5,1:5]

##            sample      TCGA_id type subtype   OR4F5
## 1 TCGA-OR-A5J1-01 TCGA-OR-A5J1  ACC      01 -9.9658
## 2 TCGA-OR-A5J2-01 TCGA-OR-A5J2  ACC      01 -9.9658
## 3 TCGA-OR-A5J3-01 TCGA-OR-A5J3  ACC      01 -9.9658
## 4 TCGA-OR-A5J6-01 TCGA-OR-A5J6  ACC      01 -9.9658
## 5 TCGA-OR-A5J7-01 TCGA-OR-A5J7  ACC      01 -9.9658

dim(mRNA_exprSet)

## [1]  7801 19729

运行EPIC

Run EPIC

#epic <- deconvo_tme(eset = expr, 
#                    method = "epic", 
#                    arrays = FALSE)
# save(epic,file = "epic.Rda")

load("epic.Rda")
dim(epic)

## [1] 7801    9

#epic[1:5,]
epic$ID <- gsub("\\.","-",epic$ID)
data <- merge(epic, 
              mRNA_exprSet, 
              by.x = "ID", by.y = "sample")
#data[1:5,]
dim(data)

## [1]  7800 19737

计算相关性

Calculate relevance

# 自定义计算相关性的函数
# Customize the function to calculate the relevance
corr <- function(gene){
  y <- as.numeric(data[, gene])
  colnames <- colnames(epic[, 2:9])
  do.call(rbind,lapply(colnames, function(x){
    dd <- cor.test(as.numeric(data[, x]), y, type = "spearman")
    data.frame(gene = gene, immune_cells = x,
               cor = dd$estimate, p.value = dd$p.value )
  }))
}
# 计算相关性
# Calculate relevance
corr("CD96")

##      gene     immune_cells          cor       p.value
## cor  CD96      Bcells_EPIC  0.201971593  1.337815e-72
## cor1 CD96        CAFs_EPIC  0.009873013  3.832944e-01
## cor2 CD96  CD4_Tcells_EPIC  0.480019990  0.000000e+00
## cor3 CD96  CD8_Tcells_EPIC  0.407591892 5.076429e-310
## cor4 CD96 Endothelial_EPIC -0.042554096  1.703373e-04
## cor5 CD96 Macrophages_EPIC  0.008188060  4.696516e-01
## cor6 CD96     NKcells_EPIC  0.067188593  2.855490e-09
## cor7 CD96  otherCells_EPIC -0.553441917  0.000000e+00

gene <- "CD96"
cancertype <- unique(data$type)
res <- data.frame(type=NULL,
                  gene=NULL,
                  immune_cells=NULL,
                  cor=NULL,
                  p.value=NULL)

for (i in cancertype){
  data2 <- data %>% filter(type==i)
  y <- as.numeric(data2[,gene])
  colnames <- colnames(epic[,2:9])
  dd2 <- do.call(rbind,lapply(colnames, function(x){
    dd <- cor.test(as.numeric(data2[,x]),y,type = "spearman")
    data.frame(type = i,
               gene = gene,
               immune_cells = x,
               cor = dd$estimate,
               p.value = dd$p.value )
             }))
  res <- rbind(res,dd2)
}
res$pvalue <- ifelse(res$p.value >= 0.05,"p≥0.05","p<0.05")
head(res)

##      type gene     immune_cells         cor      p.value pvalue
## cor   GBM CD96      Bcells_EPIC  0.42256581 1.619400e-07 p<0.05
## cor1  GBM CD96        CAFs_EPIC -0.04028372 6.340841e-01 p≥0.05
## cor2  GBM CD96  CD4_Tcells_EPIC -0.17117552 4.166795e-02 p<0.05
## cor3  GBM CD96  CD8_Tcells_EPIC  0.34204163 3.096633e-05 p<0.05
## cor4  GBM CD96 Endothelial_EPIC -0.25386741 2.299519e-03 p<0.05
## cor5  GBM CD96 Macrophages_EPIC  0.15633575 6.317888e-02 p≥0.05

write.csv(res,"output_res.EPIC.csv")

画图

EPIC只有8种细胞，画在一起就可以了。

保持跟上面TIMER组合图的风格一致。

Drawing

EPIC has only 8 types of cells, and it is enough to draw them together.

Keep the same style as the TIMER combo chart above.

ggplot(res, aes(immune_cells,
               type,
               shape = pvalue,
               color = cor)) +  
geom_point(size = 3) + 
scale_shape_manual(values = c(15,7)) + 
scale_color_gradient2(low = "#2b8cbe",
                      mid = "white",
                     high = "#e41a1c",
                   limits = c(-1, 1)) + 
theme_bw() + 
theme(axis.title.y = element_blank(),
      axis.title.x = element_blank(),
      axis.text.x = element_text(angle = 90,
                                 hjust = 1,
                                 vjust = 0,
                                family = "arial"))

ggsave("EPIC.pdf", width = 4, height = 7)

Cibersort

输入数据同“EPIC”的输入数据

运行Cibersort

Cibersort

The input data is the same as the input data of “EPIC”.

Run Cibersort

# 用IOBR包计算
# Calculated using IOBR packets
#cibersort <- deconvo_tme(eset = expr, 
#                      method = "cibersort", 
#                      arrays = FALSE, 
#                      perm = 200 )

# pancancer的免疫细胞浸润信息可以用IOBR计算，也可以直接导入官方数据 [CIBERSORT immune fractions](https://gdc.cancer.gov/about-data/publications/panimmune)
# 这里直接导入官方数据
# The immune cell infiltration information of pancancer can be calculated using IOBR or directly imported into the official data [CIBERSORT immune fractions](https://gdc.cancer.gov/about-data/publications/panimmune)
# Import official data directly here
cibersort <- read.table("TCGA.Kallisto.fullIDs.cibersort.relative.tsv", header = T)

cibersort$SampleID <- substring(cibersort$SampleID,1,15)
cibersort$SampleID <- gsub("\\.","-",cibersort$SampleID)

cibersort <- merge(cibersort,
                   mRNA_exprSet,
                   by.x = "SampleID",
                   by.y = "sample")

计算相关性

Calculate relevance

gene <- "CD96"
res <- data.frame(type=NULL,
                  gene=NULL,
                  immune_cells=NULL,
                  cor=NULL,
                  p.value=NULL)

for (i in cancertype){
  data2 <- cibersort %>% filter(type==i)
  y <- as.numeric(data2[,gene])
  colnames <- colnames(cibersort[,3:24])
  dd2 <- do.call(rbind,lapply(colnames, function(x){
    dd <- cor.test(as.numeric(data2[,x]),y,type = "spearman")
    data.frame(type = i,
               gene = gene,
               immune_cells = x,
               cor = dd$estimate,
               p.value = dd$p.value )
  }))
  res <- rbind(res,dd2)
}

## Warning in cor(x, y): the standard deviation is zero
## Warning in cor(x, y): the standard deviation is zero
## Warning in cor(x, y): the standard deviation is zero
## Warning in cor(x, y): the standard deviation is zero
## Warning in cor(x, y): the standard deviation is zero
## Warning in cor(x, y): the standard deviation is zero
## Warning in cor(x, y): the standard deviation is zero
## Warning in cor(x, y): the standard deviation is zero

res$pvalue <- ifelse(res$p.value >= 0.05,"p≥0.05","p<0.05")
head(res)

##      type gene               immune_cells         cor    p.value pvalue
## cor   GBM CD96              B.cells.naive  0.03270965 0.69814824 p≥0.05
## cor1  GBM CD96             B.cells.memory -0.02616166 0.75644282 p≥0.05
## cor2  GBM CD96               Plasma.cells  0.07045307 0.40307529 p≥0.05
## cor3  GBM CD96                T.cells.CD8 -0.18280972 0.02886319 p<0.05
## cor4  GBM CD96          T.cells.CD4.naive  0.02621328 0.75597738 p≥0.05
## cor5  GBM CD96 T.cells.CD4.memory.resting  0.06235925 0.45937045 p≥0.05

write.csv(res, "output_res.Cibersort.csv")

画图

Cibersort有22种细胞，也可以画在一起。并保持跟以上两种方法的画图风格一致。

Drawing

Cibersort has 22 types of cells, which can also be drawn together. And keep the drawing style consistent with the above two methods.

ggplot(res, aes(immune_cells,
            type,
            shape = pvalue,
            color = cor)) +  
  geom_point(size = 3) + 
  scale_shape_manual(values = c(15,7)) + 
  scale_color_gradient2(limits = c(-1, 1),
                        low = "#2b8cbe",
                        mid = "white",
                        high = "#e41a1c") + 
  theme_bw() + 
  theme(axis.title.y= element_blank(),
        axis.title.x = element_blank(),
        axis.text.x = element_text(angle = 90,
                                   hjust = 1,
                                   vjust = 0,
                                   family = "arial"))

## Warning: Removed 8 rows containing missing values or values outside the scale range
## (`geom_point()`).

ggsave("Cibersort.pdf", width = 5, height = 7)

## Warning: Removed 8 rows containing missing values or values outside the scale range
## (`geom_point()`).

Session Info

sessionInfo()

## R version 4.4.1 (2024-06-14 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 22631)
## 
## Matrix products: default
## 
## 
## locale:
## [1] LC_COLLATE=Chinese (Simplified)_China.utf8 
## [2] LC_CTYPE=Chinese (Simplified)_China.utf8   
## [3] LC_MONETARY=Chinese (Simplified)_China.utf8
## [4] LC_NUMERIC=C                               
## [5] LC_TIME=Chinese (Simplified)_China.utf8    
## 
## time zone: Asia/Shanghai
## tzcode source: internal
## 
## attached base packages:
## [1] grid      stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] IOBR_0.99.0            GSVA_1.52.3            e1071_1.7-16          
##  [4] preprocessCore_1.66.0  limSolve_1.5.7.1       limma_3.60.4          
##  [7] lubridate_1.9.3        forcats_1.0.0          stringr_1.5.1         
## [10] purrr_1.0.2            readr_2.1.5            tidyr_1.3.1           
## [13] tidyverse_2.0.0        survminer_0.5.0        clusterProfiler_4.12.6
## [16] tidyHeatmap_1.8.1      ComplexHeatmap_2.20.0  survival_3.6-4        
## [19] dplyr_1.1.4            tibble_3.2.1           EPIC_1.1.7            
## [22] showtext_0.9-7         showtextdb_3.0         sysfonts_0.8.9        
## [25] patchwork_1.3.0        ggpubr_0.6.0           ggplot2_3.5.2         
## 
## loaded via a namespace (and not attached):
##   [1] SpatialExperiment_1.16.0    IRanges_2.38.1             
##   [3] R.methodsS3_1.8.2           GSEABase_1.66.0            
##   [5] goftest_1.2-3               HDF5Array_1.32.1           
##   [7] Biostrings_2.72.1           vctrs_0.6.5                
##   [9] spatstat.random_3.3-2       digest_0.6.37              
##  [11] png_0.1-8                   shape_1.4.6.1              
##  [13] proxy_0.4-27                ggrepel_0.9.6              
##  [15] deldir_2.0-4                parallelly_1.39.0          
##  [17] magick_2.8.5                MASS_7.3-60.2              
##  [19] reshape2_1.4.4              httpuv_1.6.15              
##  [21] foreach_1.5.2               BiocGenerics_0.50.0        
##  [23] qvalue_2.36.0               withr_3.0.2                
##  [25] xfun_0.53                   ggfun_0.1.7                
##  [27] memoise_2.0.1               gson_0.1.0                 
##  [29] systemfonts_1.1.0           ragg_1.3.3                 
##  [31] tidytree_0.4.6              zoo_1.8-12                 
##  [33] GlobalOptions_0.1.2         pbapply_1.7-2              
##  [35] R.oo_1.27.0                 Formula_1.2-5              
##  [37] KEGGREST_1.44.1             promises_1.3.0             
##  [39] httr_1.4.7                  rstatix_0.7.2              
##  [41] rhdf5filters_1.16.0         globals_0.16.3             
##  [43] fitdistrplus_1.2-1          rhdf5_2.48.0               
##  [45] rstudioapi_0.17.1           UCSC.utils_1.0.0           
##  [47] miniUI_0.1.1.1              generics_0.1.3             
##  [49] DOSE_3.30.5                 S4Vectors_0.42.1           
##  [51] zlibbioc_1.50.0             ScaledMatrix_1.12.0        
##  [53] ggraph_2.2.1                polyclip_1.10-7            
##  [55] GenomeInfoDbData_1.2.12     quadprog_1.5-8             
##  [57] SparseArray_1.4.8           xtable_1.8-4               
##  [59] doParallel_1.0.17           evaluate_1.0.1             
##  [61] S4Arrays_1.4.1              hms_1.1.3                  
##  [63] glmnet_4.1-8                GenomicRanges_1.56.1       
##  [65] irlba_2.3.5.1               colorspace_2.1-1           
##  [67] ROCR_1.0-11                 reticulate_1.39.0          
##  [69] spatstat.data_3.1-2         magrittr_2.0.3             
##  [71] lmtest_0.9-40               later_1.3.2                
##  [73] viridis_0.6.5               ggtree_3.12.0              
##  [75] lattice_0.22-6              spatstat.geom_3.3-3        
##  [77] future.apply_1.11.3         XML_3.99-0.17              
##  [79] scattermore_1.2             shadowtext_0.1.4           
##  [81] cowplot_1.1.3               matrixStats_1.4.1          
##  [83] RcppAnnoy_0.0.22            class_7.3-22               
##  [85] pillar_1.9.0                nlme_3.1-164               
##  [87] iterators_1.0.14            beachmat_2.20.0            
##  [89] compiler_4.4.1              RSpectra_0.16-2            
##  [91] stringi_1.8.4               tensor_1.5                 
##  [93] SummarizedExperiment_1.34.0 dendextend_1.18.1          
##  [95] plyr_1.8.9                  crayon_1.5.3               
##  [97] abind_1.4-8                 gridGraphics_0.5-1         
##  [99] locfit_1.5-9.10             sp_2.1-4                   
## [101] graphlayouts_1.2.0          bit_4.5.0                  
## [103] fastmatch_1.1-4             textshaping_0.4.1          
## [105] codetools_0.2-20            BiocSingular_1.20.0        
## [107] pec_2023.04.12              bslib_0.8.0                
## [109] GetoptLong_1.0.5            plotly_4.10.4              
## [111] mime_0.12                   splines_4.4.1              
## [113] circlize_0.4.16             Rcpp_1.0.13-1              
## [115] fastDummies_1.7.4           sparseMatrixStats_1.16.0   
## [117] knitr_1.49                  blob_1.2.4                 
## [119] utf8_1.2.4                  clue_0.3-65                
## [121] fs_1.6.5                    listenv_0.9.1              
## [123] ggsignif_0.6.4              ggplotify_0.1.2            
## [125] Matrix_1.7-0                statmod_1.5.0              
## [127] tzdb_0.4.0                  lpSolve_5.6.22             
## [129] tweenr_2.0.3                pkgconfig_2.0.3            
## [131] tools_4.4.1                 cachem_1.1.0               
## [133] RSQLite_2.3.7               viridisLite_0.4.2          
## [135] DBI_1.2.3                   numDeriv_2016.8-1.1        
## [137] fastmap_1.2.0               rmarkdown_2.29             
## [139] scales_1.3.0                ica_1.0-3                  
## [141] Seurat_5.1.0                broom_1.0.7                
## [143] sass_0.4.9                  dotCall64_1.2              
## [145] graph_1.82.0                carData_3.0-5              
## [147] RANN_2.6.2                  farver_2.1.2               
## [149] tidygraph_1.3.1             scatterpie_0.2.4           
## [151] yaml_2.3.10                 MatrixGenerics_1.16.0      
## [153] cli_3.6.3                   stats4_4.4.1               
## [155] leiden_0.4.3.1              lifecycle_1.0.4            
## [157] uwot_0.2.2                  Biobase_2.64.0             
## [159] mvtnorm_1.3-2               lava_1.8.0                 
## [161] backports_1.5.0             annotate_1.82.0            
## [163] BiocParallel_1.38.0         timechange_0.3.0           
## [165] gtable_0.3.6                rjson_0.2.23               
## [167] ggridges_0.5.6              progressr_0.15.0           
## [169] parallel_4.4.1              ape_5.8                    
## [171] jsonlite_1.8.9              RcppHNSW_0.6.0             
## [173] bit64_4.5.2                 Rtsne_0.17                 
## [175] yulab.utils_0.1.8           spatstat.utils_3.1-1       
## [177] SeuratObject_5.0.2          jquerylib_0.1.4            
## [179] GOSemSim_2.30.2             survMisc_0.5.6             
## [181] spatstat.univar_3.1-1       R.utils_2.12.3             
## [183] lazyeval_0.2.2              shiny_1.9.1                
## [185] htmltools_0.5.8.1           enrichplot_1.24.4          
## [187] KMsurv_0.1-5                GO.db_3.19.1               
## [189] sctransform_0.4.1           rappdirs_0.3.3             
## [191] glue_1.8.0                  timereg_2.0.6              
## [193] spam_2.11-0                 httr2_1.0.6                
## [195] XVector_0.44.0              treeio_1.28.0              
## [197] gridExtra_2.3               igraph_2.1.2               
## [199] R6_2.5.1                    SingleCellExperiment_1.26.0
## [201] DESeq2_1.44.0               labeling_0.4.3             
## [203] km.ci_0.5-6                 cluster_2.1.6              
## [205] Rhdf5lib_1.26.0             aplot_0.2.3                
## [207] GenomeInfoDb_1.40.1         DelayedArray_0.30.1        
## [209] tidyselect_1.2.1            ggforce_0.4.2              
## [211] car_3.1-3                   AnnotationDbi_1.66.0       
## [213] future_1.34.0               rsvd_1.0.5                 
## [215] munsell_0.5.1               KernSmooth_2.23-24         
## [217] data.table_1.16.4           htmlwidgets_1.6.4          
## [219] fgsea_1.30.0                RColorBrewer_1.1-3         
## [221] rlang_1.1.4                 spatstat.sparse_3.1-0      
## [223] spatstat.explore_3.3-3      fansi_1.0.6                
## [225] timeROC_0.4                 prodlim_2024.06.25

FigureYa234panImmune

Houshi Xu

2025-5-20

Academic Citation

需求描述

Requirement Description

Application Scenarios

环境设置

Environment settings

TIMER

从TIMER下载数据

TIMER

Download data from TIMER

批量读入TIMER结果

Batch reading TIMER results

画图前的数据预处理

Pre-processing of data before drawing

画图

Drawing

EPIC

EPIC

输入数据

Enter data

运行EPIC

Run EPIC

计算相关性

Calculate relevance

画图

Drawing

Cibersort

运行Cibersort

Cibersort

Run Cibersort

计算相关性

Calculate relevance

画图

Drawing

Session Info